import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d
plt.rcParams['font.sans-serif'] = ['SimHei']  # 用黑体显示中文
plt.rcParams['axes.unicode_minus'] = False  # 解决负号 '-' 显示为方块的问题

def read_truth_data(filename):
    data = np.loadtxt(filename)
    wavelengths = data[:, 0]
    intensities = data[:, 1]
    return wavelengths, intensities

def read_reconstructed_csv(filename):
    # 如果csv有表头，默认用第一列和第二列
    df = pd.read_csv(filename)
    # 假设第一列是波长，第二列是强度，如果列名不确定，可以用 df.columns 查看
    wavelengths = df.iloc[:, 0].values
    intensities = df.iloc[:, 1].values
    return wavelengths, intensities

def main():
    truth_file = 'truth_spectrum.txt'
    recon_file = 'reconstructed.csv'

    truth_wl, truth_intensity = read_truth_data(truth_file)
    recon_wl, recon_intensity = read_reconstructed_csv(recon_file)

    # 对真值数据插值到重建光谱波长点
    interp_func = interp1d(truth_wl, truth_intensity, kind='linear', fill_value='extrapolate')
    truth_interp = interp_func(recon_wl)

    # 计算均方误差
    mse = np.mean((truth_interp - recon_intensity) ** 2)
    print(f"均方误差(MSE): {mse:.3e}")

    # 绘制对比图
    plt.figure(figsize=(8, 5))
    plt.plot(recon_wl, recon_intensity, label='重建光谱', linewidth=2)
    plt.plot(recon_wl, truth_interp, linestyle='--', label='真值光谱（插值后）', linewidth=2)
    plt.scatter(truth_wl, truth_intensity, color='red', s=20, label='真值光谱原始点')
    plt.xlabel("波长 (nm)")
    plt.ylabel("强度 (a.u.)")
    plt.title("重建光谱与真值光谱对比")
    plt.legend()
    plt.grid(True)
    plt.tight_layout()
    plt.show()

if __name__ == '__main__':
    main()
